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Abstract 

In  this  paper  a  Boundary  Integral  Equation  Method  (BIEM)  is  described 
for  the  computation  of  scattering  from  a  finite,  rigid,  cylinder  near  a  pressure- 
release  interface.  The  cylinder  lies  parallel  or  tilted  with  respect  to  the  interface 
plane  so  that  the  azimuthal  symmetry  of  the  problem  is  destroyed.  The  scatter¬ 
ing  solution  is  first  described  in  terms  of  an  azimuthally-symmetric  free  space 
solution.  The  multiple  interactions  of  the  scattered  field  with  the  interface  are 
accounted  for  by  an  azimuthal-conversion  matrix.  In  the  numerical  examples, 
the  method  of  this  paper  is  benchmarked  using  wavefield  superposition  for  a 
sphere  near  a  pressure-release  surface.  Computed  scattered  spectra  are  shown 
for  a  finite  cylinder,  parallel  and  tilted  with  respect  to  the  interface,  and  for 
a  variety  of  source/receiver  geometries.  The  differences  resulting  from  not  in¬ 
cluding  multiple  target/interface  interactions  (single  scatter  solution)  and  from 
including  all  interactions  are  presented.  The  problem  of  irregular  frequencies 
for  the  single-scatter  and  the  fully  coupled  BIEM  are  discussed  and  numerically 
examined. 


PACS  numbers:  43.20.Fn,  43.30.Gv 
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I.  INTRODUCTION 

The  numerical  solution  of  scattering  from  objects  near  or  on  an  interface  has  been 
described1^6  by  previous  authors.  In  Refs.  4  and  5,  the  scattering  was  solved  using 
a  finite-element  code.  For  some  geometries,  such  as  a  sphere  or  a  vertically-oriented 
cylinder  near  an  interface,  the  scattering  problem  can  be  solved  as  a  sequence  of  az¬ 
imuthal  problems  (i.e. ,  the  Fourier  components  are  uncoupled).  For  a  cylinder  lying 
parallel  or  tilted  with  respect  to  the  interface,  the  azimuthal  symmetry  is  broken  and 
one  utilizes  either  a  single-scattering  solution  or  a  full  three-dimensional  solution5. 
The  single-scattering  solution  is  an  approximate  solution  where  the  effects  of  the 
interface  are  accounted  for  with  respect  to  the  incident  field  and  also  in  the  propa¬ 
gation  of  the  field,  computed  on  the  surface  of  the  cylinder,  to  a  receiver (s).  For  the 
single-scatter  approach  the  solution  upon  the  target’s  surface  can  still  be  computed 
with  an  uncoupled  azimuthal  expansion  method.  The  accuracy  of  the  single  scatter 
solution  has  been  explored  by  various  authors2,6  and  although,  it  is  often  a  very  good 
approximation,  it  is  also  inaccurate  for  certain  scenarios.  This  paper  will  consider 
the  accuracy  of  the  single-scatter  approximation  for  different  monostatic  scattering 
angles  and  different  orientations  of  a  rigid  target  near  a  pressure-release  surface.  A 
flat-ended  cylinder,  tilted  or  parallel  to  the  interface,  serves  as  an  interesting  test 
case  for  the  single-scatter  approximation.  The  scattering  will  consist  of  both  spec¬ 
ular  and  diffraction  returns.  As  well,  for  a  tilted  cylinder,  there  can  be  significant 
target /interface  interactions  for  certain  geometries. 

The  solution  method,  Boundary  Integral  Equation  Method  (BIEM),  using  Fourier 
expansions  is  first  described  for  the  free  space  problem7,8.  Next,  a  finite,  rigid  cylinder 
below  a  pressure  release  surface  is  considered.  In  this  case,  the  concept  of  an  image 
cylinder  is  introduced  and  a  conversion  matrix,  describing  the  conversion  of  the  so¬ 
lution  upon  this  image  cylinder  into  a  field  re-incident  upon  the  physical  cylinder  is 
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derived.  This  conversion  matrix  can  be  ignored  in  the  BIEM  system  of  equations 
(leading  to  the  single-scatter  solution),  can  be  included  iteratively,  or  a  fully-coupled 


system  of  equations  (all  orders  of  interactions)  can  be  solved.  For  the  cylinder  ly¬ 


ing  parallel  to  the  interface,  a  computationally  efficient  method  for  computing  the 
conversion  matrix  is  derived.  For  a  cylinder  tilted  with  respect  to  the  interface,  the 
conversion  matrix  is  somewhat  more  difficult  to  compute. 

In  the  Numerical  Examples  section,  the  case  of  a  sphere  lying  below  a  reflecting 
interface  is  first  considered  as  a  benchmark  case.  With  respect  to  the  vertical  axis,  the 
problem  is  azimuthally  symmetric  and  we  solve  it  with  a  virtual  source  method9-11 
using  the  appropriate  half-space  Green’s  function.  However,  we  also  consider  the 
Fourier  expansion  about  the  horizontal  axis  and  use  the  cylindrical  BIEM  with  a 
conversion  matrix.  Next,  using  the  BIEM  technique,  a  finite  cylinder  lying  parallel 
to  the  interface  is  considered.  The  spectra  for  the  scattered  fields  are  computed 
for  3  monstatic  angles  using  the  single-scatter  and  full-interaction  solutions,  and  the 
differences  are  discussed.  The  irregular  frequencies7,8,12  of  the  single-scatter  and  full 
multiple-scatter  solutions  are  investigated.  Finally,  a  cylinder  tilted  at  -45°  with 
respect  to  the  horizontal  is  considered.  Once  again,  the  single-  and  full-scattering 
solutions  are  computed  for  3  monstatic  angles  and  the  spectra  and  corresponding 
time  series  computations  discussed. 

II.  THEORY 

A.  A  cylindrical  target  in  free  space 

The  BIEM  for  a  cylinder  in  free  space  is  first  described.  First,  the  problem  is  con¬ 
sidered  at  a  fixed  angular  frequency  c o  and  we  write  the  three-dimensional  Helmholtz 
Equation,  for  a  rigid  scatterer,  as 


(1) 
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Figure  1:  Schematic  drawing  of  a  1.5m  long  cylinder  with  radius  0.5.  A  generic 
receiver  position  is  shown  at  (r,  6,z  =  0).  A  point  on  the  surface  at  (r1  =  0.5,  6',  z’  =  0) 
is  also  shown. 


Here  we  are  taking  the  normal  vector  to  point  out  of  the  cylinder  into  the  surrounding 
fluid.  The  coordinates  r,  r'  are  the  radial  coordinates  with  respect  to  the  cylinder  axis 
and  z,  z!  are  the  coordinates  along  the  axis  of  the  cylinder.  The  normal  vector  will 
have  components  in  the  r’  and  zl  directions  of  the  cylinder.  A  schematic  of  a  three- 
dimensional  cylinder  and  its  notation  is  shown  in  Fig.  1. 

The  three-dimensional  Green’s  function  G  is  given  by 


G(r,z,F,r',z',F)  = 


exp(iuj/cjR  —  R  '|) 
Att\R-R  '| 


(2) 
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where  c  is  the  speed  of  sound  in  the  surrounding  fluid  and  R  denotes  the  three- 
dimensional  position  vector.  Let  us  expand  the  unknown  pressure  hied  on  the  surface 
of  the  target,  p(r,  z,  9),  as 


M 

p(r ,z,9)  —  Pm(r>  z)  exp(im0)  (3) 

m=—M 

where  the  infinite  Fourier  sum  has  been  truncated. 

We  write  the  distance  term  in  the  Green’s  function  of  Eq.(2)  as 


R  —  R  ' |  =  \Jr 2  +  (r')2  —  2 rr'cos{9  —  9’)  +  (z  —  z')2. 


(4) 


Substituting  the  Fourier  expansion  for  p,  Eq.(3)  into  Eq.(l)  ,  the  expression 

r  r2n  q  exp(hn/m/r2  +  (r')2  —  2 rr'cos(9  —  O’)  +  (z  —  z')2) 

p(r,z,0)  =  pmc+  /  /  — 7 - /  = - 

JS'  Jo  on  An yjr2  +  (r')2  —  2 rr'cos{9  —  9')  +  (z  —  z')2 


M 


x  ^2  Pm{r\  z')  exp{im6')dO'r\S')dS' 

m=—M 


(5) 


is  obtained.  The  (r',  z')  perimeter  of  the  cylinder  is  denoted  by  S'.  Performing  the 
azimuthal  integral  and  using  the  substitution  that  ip  =  (9  —  6')  one  obtains 

M 

p(r,z,0)  =  pinc  +  exp (im0)  gm(r,  z\r' ,  z')pm{r' ,  z')r' (S')dS'  (6) 

rnt^M  JS' 


where 


r2ir  Q  exp(hn / cJ r2  +  (r')2  —  2rr'  cos  (ip)  +  (z  —  z')2) 
gm(r,z;r',z')=  j  — - ; - ,  [  . . ^ - exp (-irmp)dip. 


I  o  dn'  4nsjr2  +  (r')2  —  2  rr'cos{ip)  +  (z  —  z')2 


(7) 


Multiplying  Eq.(6)  by  exp(— imi9)/(2n)  and  integrating  with  respect  to  9  over  the 
interval  [0  2i r],  we  obtain  from  the  orthogonality  of  exp(—  irrii9) 


Pmi(r,z)  =p™ci(r,z )  +  /  (8) 

J  S' 


DRDC-RDDC-20 1 4-P26 


Fawcett-  J AS A 


7 


where 

Pm t(r,  z)  =  ^  Jo  exp(-2mi 6)pmc(r,  z,  6)dO.  (9) 

(it  should  be  noted  that,  in  fact,  we  will  compute  both  g_mi  and  gmi  in  Eq.(7)  using 
a  cosine  transform  over  [0, 7r]).  Thus,  for  each  angular  order  m\  we  have  a  system  of 
equations  to  solve, 

=»Cr  (10) 

where  Rmipmi  represents  a  discrete-matrix  version  of 


'S' 


5lm1(r,  Z'i  r'(S'),  z'(S'))pmi(r\S'),  z'(S'))r'(S')dS' . 


(11) 


The  matrix  /  in  Eq.(10)  is  the  identity  matrix  and  the  factor  1/2  is  a  result  of  the 
standard  limit  for  the  Helmholtz  equation  as  the  observation  point  approaches  the 
surface. 

In  order  to  discretize  the  integral  of  Eq.(ll)  we  discretize  the  contour  S'  into 
discrete  panels.  We  consider  pmi  to  be  constant  over  a  panel  and  represent  the  value 
at  the  midpoint  of  the  panel.  Similarly  we  use  the  midpoint  value  of  r' .  For  gmi 
we  will  consider  the  midpanel  values  of  (r,  z)  and  (r',  z')  except  for  cases  where  gmi 
exhibits  singular  behaviour  and  we  will  perform  an  analytic  integration  over  the  panel. 

From  Eq.(7)  ,  there  is  a  singular  behaviour  to  gm  as  (r',  z')  — >  (r,  z).  Let  us  rewrite 
Eq.(7)  as 


f  f  r^-rr  Q  r  exp(hn/cwr2  +  (r')2  —  2rr'cos{^)  +  [z  —  z')2) 

gm(r,z;r,z )  =  /  —  < - ,  s  - exp(im^) 

Jo  on  l  47t yr2  +  (r')2  —  2rr,cos('ijj)  +  (z  —  z’)2 

1  ^^+Gs(r,z-  r’,z’)  (J2) 


Air^Jr2  +  (r')2  —  2  rr'cosfy)  +  (z  —  z')2 


4n 


where 


Gs(r,z,;r', 


_d _ 

dn'  ry, 2  _|_  (jr'^  2 


2 rr'cos(^)  +  (z  —  z')2  ' 


(13) 
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The  integral  of  Eq.(13)  can  be  expressed  analytically11.  In  particular,  following 
the  notation  of  Ref.  11, 


Gs(r,  z,  \r',z')  =  2 nr> 


{R2[E{k)  -  AT(«)(1  -  k2)]  -  2 r'(r'  +  r)E(n)} 
r'RK  1  —  k2) 


4(z'  —  z)E(k) 

Wi  -  *2) 


(14) 


where  nr>  and  nz>  are  the  2  normal  components  along  the  cylinder  profile  and  K  and 
E  are  the  complete  elliptic  integrals  of  the  first  and  second  kinds  respectively.  The 
terms  Re  and  k  are  defined  by 


(15) 


The  first  term  in  Eq.(14)  is  singular  as  r  r'  ,  z  z' ,  and  k  — >  1.  However,  this 
term  can  be  analytically  integrated  over  the  panel  where  k  =  1  by  using  an  expansion 
of  K(k)  in  the  neighbourhood  of  k  =  1. 

In  the  case  of  a  flat-ended  cylinder,  there  is  another  source  of  singular- like  terms 
which  occur  in  the  matrix  kernel  for  points  near  the  edges  of  the  flat  endcaps.  Con¬ 
sider,  for  example,  the  second  term  in  Eq.(14)  for  ( r',z ')  on  the  flat  endcap  and 
(r,  z)  along  the  length  of  the  cylinder.  On  the  endcap  nz/  =  1,  but  since  z  is  along 
the  length  of  the  cylinder,  z'  —  z  =  Az  is  non-zero  with  the  denominator  approach¬ 
ing  zero  as  (■ r',z ')  —>■  ( r,  z ).  For  these  cases,  (we  will  consider  the  discrete  value  of 
r1  =  rj )  one  can  subtract  a  term  with  the  correct  singular  behaviour  and  add  back  in 
a  corresponding  analytical  integral, 


r  -  rj  -  5r/ 2 


z  —  z' 


))  (16) 
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Similarly,  the  nr>  term  in  Eq.(14)  gives  rise  to  an  almost  singular  term  when  (■ r',z ') 
is  on  the  length  of  the  cylinder  and  (r,  z)  is  on  an  endcap.  In  order  to  deal  with 
these  behaviours  near  the  endcaps,  simple  expressions  which  capture  the  singular 
behaviours  are  subtracted  from  all  the  terms  of  Eq.(14)  and  then  analytical  integrals 
of  these  terms  are  added  back  in. 

The  BIEM  approach  of  this  paper  suffers  from  the  problem  of  irregular  frequencies7,8,12 
that  is,  specific  frequencies  where  the  BIEM  system  of  equations  does  not  have  an 
unique  solution.  In  the  free  space  case,  these  frequencies  correspond  to  the  inte¬ 
gral  equation  (the  continuous  version  of  Eq.(10))  having  a  zero  eigenvalue.  This,  in 
turn,  corresponds  to  a  zero  eigenvalue  for  the  Dirichlet  problem  (i.e.  p  =  0  on  the 
boundaries)12  for  the  interior  of  the  target.  For  the  flat-ended  cylinder,  these  irreg¬ 
ular  frequencies  can  be  analytically  predicted.  The  interior  Dirichlet  problem,  for 
azimuthal  order  m,  is  solved  by  eigenfunctions  of  the  form 

T(r,  z)  =  sin(nn  /  Lz)  Jm((y4a2  /  c2  —  ( nir/L)2)r ).  (17) 

Denoting  the  radius  of  the  cylinder  as  ro,  then  for  a  given  azimuthal  order,  the 
irregular  frequencies  are  given  by 

fn,g  =  +  (“)2-  rc=l,...,oo,g  =  l,...oo  (18) 

(see,  also,  Ref.  12)  where  rq  is  a  zero  of  Jm{u).  In  this  paper,  we  will  add  Nsh  extra 
interior  points  at  which  to  evaluate  the  Helmholtz  equation.  As  these  field  points  are 
in  the  interior  of  the  cylindrical  object,  the  interior  Helmholtz  equation  applies.  For 
example,  for  the  set  of  interior  points  at  (rq,zq),  the  interior  equation  is: 

N 

ZC(r9>  Zq)  =  -J2  Rm(rq,  Zq,  Tj ,  Zj)p(j)  (19) 

3= 1 

where  Rrn  is  the  discrete  kernel  at  angular  order  m.  In  our  implementation,  the  Nsh 
points  are  displaced  along  radials  connecting  the  interior  origin  of  the  cylinder  to  its 
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(r',  z')  perimeter.  The  distance  to  be  displaced  along  the  radials  is  an  user  input. 
For  each  radial,  this  distance  is  varied  by  a  small  random  amount.  For  each  angular 
order,  m,  we  will  have  a  ( N  +  Nsh)  x  N  system  of  equations  to  solve  in  a  least-squares 
sense.  As  was  shown  in  Ref.  12,  the  resulting  system  of  equations  now  has  an  unique 
solution. 


B.  A  cylinder  below  an  interface 


The  previous  derivation  of  the  BIEM  considered  an  axi-symmetric  target  in  free 
space.  In  the  case  that  the  target  is  below  a  pressure-release  interface,  the  appropriate 
three-dimensional  Green’s  function  is 


G(R-R') 


exp(iu / c\R  —  R  1 1)  exp(iu / c\R  —  R  '|) 

in\R-R'\  iir\R-R'\ 


(20) 


where  R  '  denotes  the  position  vector  of  the  image  point  of  R' . 

This  expression  can  be  used  in  Eq.(l)  .  Once  again,  the  unknown  pressure  held  on 
the  surface  of  the  cylindrical  object  can  be  expanded  in  a  Fourier  series  with  unknown 
coefficients.  However,  in  general,  the  sets  of  equations  for  the  various  azimuthal  orders 
are  coupled,  unlike  the  free  space  case  where  the  integral  equations  for  the  different 
azimuthal  orders  were  imcoupled.  The  one  exception  occurs  when  the  cylinder  is 
oriented  vertically  so  that  the  vertical  axis  is  the  axis  of  symmetry  for  the  scattering 
problem.  In  this  case,  the  systems  of  equations  for  the  various  azimuthal  orders 
remain  uncoupled. 

Instead  of  using  the  full  three-dimensional  Green’s  function  of  Eq.(20)  in  Eq.(l) 
to  account  for  the  presence  of  the  interface,  let  us  consider  an  iterative  solution  to 
the  problem  using  the  free-space  Green’s  function.  The  zero’th  order  solution  So  will 
be  determined  by  considering  the  incident  held  from  a  point  source  and  its  negative 
image  with  respect  to  the  upper  interface.  However,  the  solution  on  the  surface  of 
the  cylinder  will  be  computed  by  using  only  the  free-space  Green’s  function  in  the 
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integral  equation.  The  azimuthally  uncoupled  approach  can  be  used  in  this  case. 
This  solution  can  then  be  propagated  back  to  a  receiver  using  the  half-space  Green’s 
function.  Thus,  the  upper  interface  is  accounted  for  with  respect  to  the  incident 
and  scattered  fields,  but  not  in  the  solution  of  the  integral  equation  itself.  This 
corresponds  to  the  single-scattering  solution  as  discussed  by  various  authors2,6. 

In  reality,  the  field  scattered  from  the  target  reflects  from  the  upper  interface  and 
is  reincident  upon  the  target.  This  rescattering  process  will  be  characterized  by  an 
azimuthal  conversion  matrix  Cnm(Rj,  Rj).  This  matrix  expresses  the  nth  azimuthal 
coefficient  of  the  pressure  field  reincident  upon  the  cylinder  at  a  discrete  point  Ri  due 
to  a  solution  coefficient  of  unity  for  the  mth  azimuthal  order  at  the  discrete  point 
Rj.  Considering  all  the  possible  azimuthal  orders  (positive  and  negative  as  well) 
and  all  the  discretization  points  of  the  cylinder,  there  is  a  global  conversion  matrix 
G  of  dimension  (2M+l)(N+Ns/j)  x  (2M+1)N  where  M  is  the  number  of  azimuthal 
orders,  N  is  the  number  of  discretization  points  and  Nsh  is  the  number  of  added 
interior  points  for  Schenke’s  method.  Let  us  denote  the  zero-til  order  solution  upon 
the  target  as  so-  Then  the  series  or  iterative  solutions  can  be  formally  expressed  as 

s  =  So  +  (Go)  1C's*o  +  (Go)  1C(Go)  lCso  +  ...((Go)  1C)Nso  +  ... 

=  (/  -  (G0)_1G)-1s0  =  (G0  -  G)_1G0s0  =  (G0  -  C)~lpinc.  (21) 

The  concept  of  irregular  frequencies,  and  their  relationship  to  the  eigenvalues  of 
the  interior  Dirichlet  problem  was  discussed  for  the  case  of  a  target  in  free  space.  In 
the  case  of  a  target  below  an  interface,  the  irregular  frequencies  for  the  exact  BIEM 
formulation  are  the  same  as  for  the  free  space  case.  In  order  to  prove  this,  the  theorems 
of  Ref.  12  (text  and  appendix)  are  followed.  However,  instead  of  considering  the  free 
space  Green’s  function  (or  its  normal  derivative)  as  the  kernel  in  these  theorems,  one 
can  use  the  half-space  Green’s  function  (with  the  target  considered  in  the  half-space 
domain)  without  changing  the  approach  of  Ref.  12. 
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C.  Computation  of  the  conversion  matrix  C 

We  now  consider  the  computation  of  the  conversion  matrix  G,  either  to  use  in 
the  iterative  approach  or  in  the  exact  matrix  (Go  —  G)  of  Eq.(21)  .  Given  a  pressure 
distribution  upon  the  object,  what  is  the  resulting  field  reincident  upon  the  object 
due  to  the  scattered  field?  We  found  it  helpful  in  our  derivation  to  replace  the  original 
half-space  problem  with  the  equivalent  infinite  space  problem,  where  there  is  negative 
image  source  reflected  above  the  interface  and  also  a  reflected  image  of  the  target 
shape.  In  the  following,  we  refer  to  the  target  in  the  water  as  the  real  cylinder  and 
the  reflected  target  as  the  image  cylinder.  There  is  a  local  (r,  z,  6)  coordinate  system 
for  each  of  the  cylinders  where  9  on  the  image  cylinder  has  the  same  orientation  as 
on  the  real  cylinder.  These  coordinate  systems  are  illustrated  in  Fig.  2.  Due  to  the 
mirror-symmetry  of  the  problem,  if  the  pressure  field  on  the  real  cylinder  is  denoted 
as  pr,  then  the  pressure  field  pj  on  the  image  cylinder  is  given  by 


=  -PR-ir'i,  4  71  -  0[) 


(22) 


where  for  pr  the  coordinates  on  the  real  cylinder  are  the  mirror-image  coordinates  of 
those  on  the  image  cylinder.  In  terms  of  the  Fourier  coefficients  cm  of  the  pressure 
field  on  the  real  cylinder, 

M 

Pi(rp,zp)=  -c-m(rp,Zp)exp(im9i)exp(imn).  (23) 

m=—M 

Now,  let  us  consider  the  field  incident  upon  the  real  cylinder  due  to  a  pressure 
field  of  the  form  of  Eq.(23)  on  the  image  cylinder.  The  three-dimensional  free-space 
Green’s  function  is  used  to  propagate  the  field  from  the  image  cylinder  back  to  the 
real  cylinder, 


Pr  =p(ri,Zi,Oi)  =  - 


r27T  dG(r i,  Zi,  Op  r',  z\,  0') 


Pide'ds' 


(24) 


I  s'  Jo  dn' 

where  all  the  coordinates  refer  to  the  coordinate  system  of  the  image  cylinder.  For 
the  representation  of  the  Green’s  function,  the  following  Fourier-Bessel  expansion  is 
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Figure  2:  Schematic  drawings  of  geometry  and  notation  for  (upper)  tilted  cylinder 
with  a  source  and  the  mirror  image.  The  tilted  cylinder  has  a  centre  offset  H.  (lower) 
The  angular  coordinate  system  is  shown  on  the  endcaps  of  a  cylinder  lying  parallel  to 
the  interface.  A  receiver  point  P  is  indicated  with  the  angles  9r  and  $i  the  azimuthal 
angles  measured  from  the  lower  and  image  cylinders,  respectively. 
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used, 


i  M 

G{ri,zi,ei\Tfi,^M  =  —  Y 

271  m=—M 


exp (ikz(zi  -  z'i))Jm('yr,i)Hm('yri) / (4i)dkz)  exp (*m(0j-0')). 


(25) 

(this  form  can  be  simply  derived  by  taking  the  azimuthal  and  z-Fonrier  transforms 
first  and  then  constructing  the  radial  Green’s  function  for  azimuthal  order  m).  For 
the  case  of  a  finite  cylinder  with  constant  radius  r0,  then  the  normal  derivative  of 
the  Green’s  function  corresponds  to  either  the  r\  derivative  along  the  length  of  the 
cylinder  and  the  z\  derivative  at  the  end  caps.  Utilizing  Eqs.(23)-(25),  and  perfoming 
the  azimuthal  integration  in  Eq.(24),  the  resulting  expression  for  the  field  reincident 
upon  the  real  cylinder  due  to  the  solution  coefficient  cm  is 


@i) 

X 


-  (^cm(r'  =  r0’z'i)  (^J_^exp(ikz(zi  -  z^J-M 
H_m(p(ri)dk^j  r0dz'^J  x  exp(— imfy)  exp(— imn)/(Ai) 


+  endcap  contributions 


(26) 


In  the  following,  we  will  concentrate  upon  the  expressions  for  the  reincident  field  due 
to  the  pressure  field  along  the  length  of  the  image  cylinder.  The  expressions  for  the 
endcap  contributions  are  similar,  but  involve  the  z-derivative  of  the  Green’s  function. 
For  example,  the  expression  for  the  contribution  from  the  endcap  (Z  =  L)  with  the 
normal  in  the  positive  z  direction  is  given  by 

rro  /  roo  \ 

Pm,ec  =  -j  cm{r'  =  r0,  z')  ij  ~ikz  exp(ikz(zi  -  L))  J-m('jr,i)H_m('yri)dkz)  ridr  ■ 
x  exp(— imOi)  exp(— zm7r)/(4z).  (27) 

In  order  to  utilize  Eqs.(26)  and  (27)  for  real  cylinder  coordinates,  they  will  need 
to  be  converted  into  the  image  cylinder  coordinates.  One  of  the  discretization  points 
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(r9,  zq )  along  the  real  cylinder  is  considered.  For  the  value  of  rq,  the  Cartesian  values 
as  a  function  of  an  angular  variable  are  considered 

xq(6)  =  rq  cos  0,  yq(9)  =  rq  sin  9  (28) 

The  Cartesian  points  (xq(9),  yq(9),  zq)  are  then  converted  into  the  coordinate  system 
of  the  image  cylinder  to  yield  the  new  coordinates  Zi(zq,xq(6)),ri(zq,xq(9),yq(9)), 
9i(zq,xq(9),yq(9)).  We  now  consider  a  discrete  version  of  Eq.(26)  for  one  of  the  the 
B1EM  discretization  points,  z\  =  zp,  r\  =  rp.  The  conversion  matrix  for  the  (n,  m) 
element  is  computed  by  projecting  exp(— inO)  upon  the  expressions  from  Eqs.(26) 
and  (27)  , 


p2i r  poo 

Cnm{rq,zq-,rp,zp)  =  -—  /  exp(—in9)  exp(—imdi(zq,  xq(9),yq(9)))  (29) 

Z7T  J 0  J — oo 

d  Tq  /\ 

x  exp (ikz(zi(zq,  xq(9))  -  zp))  —  J-m(jrp)H-rn('yri(zq,  xq(6),yq(9)))  exp(irmr)  —^-dOdkz. 

where  Ap  is  the  discrete  length  of  the  pt\i  B1EM  panel.  In  this  equation  the  term 
cm(r'  =  r0,z')  of  Eq.(26)  is  replaced  by  the  value  unity  as  we  are  considering  the 
conversion  matrix. 

The  azimuthal  integral  of  Eq.(29)  is  performed  by  utilizing  a  numerical  discretiza¬ 
tion.  The  expression  of  Eq.(29)  is  for  (rp,  zp)  on  the  length  of  the  cylinder.  There  is  a 
simlar  expression  for  (rp,  zp)  on  the  endcaps.  Numerically,  all  the  (p,  q )  sets  of  points 
are  considered  one-by-one  and  the  coefficients  for  all  values  of  the  input  and  output 
azimuthal  orders  ( n,m )  computed.  The  expression  of  Eq.(29)  is  not  always  valid  as 
it  is  assumed  in  the  derivation  that  Ty  >=  rp.  For  example,  if  the  cylinder  is  oriented 
vertically  there  are  values  of  ? y  on  the  endcaps  and  values  of  rp  on  the  image  cylinder 
such  that  7y  <  rp  and  this  representation  breaks  down.  There  are,  however,  alternate 
valid  representations  which  can  be  used  in  this  case. 
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In  the  case  that  the  real  cylinder  is  lying  parallel  to  the  interface,  the  mathematics 
can  be  simplified  using  the  addition  theorem13  for  cylindrical  Bessel  Functions.  In 
particular,  in  Eq.(26),  the  term 


OO 

H~m('yri)exp(-im9i)  =  exp(-imn)  ^  Jn(ir)F[n+m(^D)  exp(in6)  exp(imr)  (30) 

n=— oo 

where  r,  6  refer  to  coordinates  measured  from  the  axis  of  the  real  cylinder  and  D  is 
twice  the  distance  from  the  interface  (H  in  Fig.  2). 

In  this  case,  one  obtains 

Cnm{r q,  zq]  r p,  zp)  =  —  exp(m7r ) / (4i)r PAP 

/oo  d 

exp (ikz(zq  -  zp))  —  Jm(rrrp)Jn(ryrq)Hn_m(^D)dkz (31) 

This  result  is  for  (rp,  zp)  along  the  length  of  the  cylinder.  The  result  for  a  discrete 
panel  in  one  of  the  endcaps  is  similar, 

Cnm(rq,zq]rp,zp)  =  —  exp(imr)/ (4i)r  pAp  (32) 

/OO 

-ikzazexp(ihz(zq  -  zp))Jm('yrp)Jn(^rq)Hn_m(^D)dkz 

-OO 

where  az  —  1  for  the  upper  endcap  and  az  =  —  1  for  the  lower  endcap.  These  integrals 
would  appear  to  involve  a  fair  amount  of  computation.  However,  the  computations 
can  be  made  very  vectorizable.  The  kz  integral  is  truncated  and  discretized  for  a  set 
of  (kZjj,j  =  1,  ...,NZ).  A  Nj  x  Nz  matrix  is  computed  for  Nj  the  number  of  discrete 
panels  for  the  BIEM.  A  matrix  can  be  defined  for  the  points  along  the  image  cylinder, 

Cs(p,  kzj)  =  -exp(imr)cm(rp,zp)ApAkz  j 

x  rpfp(ljrp )  f^p(~'ikz,jZp)Hn_m('YjD)  dJm^rp)  (33) 

and  a  matrix  for  points  along  the  real  cylinder, 

CR{q,  kzJ)  =  exp(ikZjZq)Jn('yjrp)/(4:i).  (34) 
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Here  kzj  is  a  discrete  value  of  kz,  A *.  is  the  corresponding  differential  (computa¬ 
tionally,  we  will  add  a  small  imaginary  component  to  kz  and  thus  the  differential  is 
complex  and  in  our  implementation  varies  with  respect  to  the  real  value  of  kz).  Then 
Eq.(32)  can  be  written  compactly  as 

Cnmij'q)  Zq,  Tp,  Zp)  Cr(^Q,  kzJ^)Cg  (p,  (35) 


It  can  be  seen  that  for  different  azimuthal  orders  (n,  m),  Cr{ji)  and  Cgirn )  are  not  cou¬ 
pled  with  respect  to  the  azimuthal  orders  with  the  exception  of  the  term  Hn_m(7JD). 

III.  NUMERICAL  EXAMPLES 


As  a  benchmark  case,  a  rigid  sphere  beneath  the  interface  is  considered.  In  this 
case,  the  most  straightforward  approach  to  solving  the  problem  is  to  consider  the  axis 
of  symmetry  as  being  the  vertical  axis.  Even  with  the  upper  interface,  the  problem 
is  still  axi-symmetric  and  one  can  use  a  de-coupled  Fourier  approach.  There  are  a 
variety  of  approaches  for  solving  this  problem.  In  this  paper,  we  use  the  method  of 
virtual  sources  or  waveheld  superposition9-11.  A  set  of  301  point  sources  of  the  form 

1  ,ex.p(ikRi)  exp  (ikR2), 


An  R\ 


R2 


(36) 


is  used,  where 


R\  =  sjr2  +  (r')2  —  2 rr'  cos^  +  (z  —  z1)2) 
R2  =  \jr2  +  (r')2  —  2 rr'  cos  ij)  +  (z  +  z')2) 


(37) 


and  if)  —  9  —  9' .  The  second  term  in  Eq.(37)  accounts  for  the  upper  interface.  The 
coordinates  for  a  set  of  point  sources  (r',  z')  are  displaced  radially  inwards  from  the 
surface  of  the  sphere  by  hi  =  0.007  m.  The  point  sources  for  the  mth  azimuthal  order 
are  computed  from  the  three-dimensional  point  sources, 

/*7 r 

sm  1  —  /  S\(r,  z;  r\  z'\  ip)  cosm^d^  (38) 

Jo 
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As  with  the  standard  BIEM,  the  method  of  wavefield  superposition  suffers  from  ir¬ 
regular  frequencies10.  As  a  remedy  to  this,  we  consider  another  set  of  point  sources, 
sm  2  displaced  inwards  from  the  sphere’s  surface  by  <52  =  0.0072  m.  Then  we  consider 
the  composite  Green’s  functions  of  the  form 


sm  =  (1  +  1000z/fc)sm,i  +  (1  -  1000i/fc)smi2. 


(39) 


This  form  of  Green’s  function  has  both  monopole  and  dipole  contributions.  For 
each  angular  order  m  we  solve  a  system  of  equations  of  the  form, 

=  i  =  i . N, j  =  1 . N  (40) 

On  on 

where  am  is  the  vector  of  unknown  coefficients  and  N  =  301  in  our  case.  Once  these 
coefficients  have  been  determined,  then  the  pressure  field  at  a  receiver  location  can 
be  constructed  by  a  Fourier  superposition  of  the  fields  computed  at  the  receiver  for 
each  azimuthal  order. 

In  the  first  example,  a  rigid  sphere  of  radius  0.4  m  is  located  1  m  below  the 
upper  scattering  surface.  We  use  both  the  superposition  and  the  cylindrical  BIEM 
(but  for  a  sphere)  methods  to  compute  the  scattered  pressure  field  as  a  function  of 
frequency  for  source/receiver  pairs  located  at  a  range  of  10m  from  the  sphere  centre 
(below  the  sphere)  and  at  angles  of  0,30,  and  60  degrees  off  the  vertical  axis  (axis  is 
oriented  downwards).  We  will  compute  the  full-scattering  solutions.  For  the  virtual 
source  technique  we  used  301  sources  at  all  frequencies  -  the  number  of  azimuthal 
orders  was  computed  from  the  formula  Ng  =  round (1.5a;  x  0.4/1500.  +  6).  For  the 
cylindrical  BIEM  approach,  we  used  201  discretization  points  along  the  surface  and 
17  additional  interior  points  for  Schenck’s  method.  The  number  of  azimuthal  orders 
was  computed  from  Ng  =  round(1.5a;  x  0.4/1500.  +  5)  (or,  in  general,  for  the  value  of 
0.4  we  use  the  radius  of  the  cylinder).  It  should  be  emphasized  that  for  the  cylindrical 
BIEM,  the  axis  of  symmetry  of  the  sphere  is  taken  to  be  the  horizontal  axis  and  the 
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Figure  3:  A  comparison  of  the  computed  backscattered  spectra  using  the  virtual 
source  technique  (sold  lines)  and  the  cylindrical  BIEM  (markers)  for  source/receiver 
60°  (circles),  30°  (diamonds),  and  0°  (+)  off  the  vertical  axis  below  the  sphere.  For 
clarity,  the  symbols  are  shown  at  every  fifth  frequency  value. 

azimuthal  problem  is  fully  coupled.  In  Fig.  3,  the  spectra  of  the  scattered  holds  as 
computed  using  the  virtual  source  (solid)  and  the  cylinder  BIEM  method  (markers) 
are  shown  for  the  angles  of  0,30,  and  60  degrees  (0°  corresponds  to  normal  incidence 
from  below).  As  can  be  seen  the  agreement  between  the  2  methods  is  excellent. 

Having  obtained  excellent  agreement  of  our  cylindrical  method  with  a  virtual 
source  approach  for  a  sphere  near  an  interface,  the  case  of  a  hat-ended  cylinder  is 
considered.  The  cylinder  is  discretized  with  301  points  and  21  interior  points  are 
used  for  Schenck’s  method.  It  has  a  total  length  of  1.5  m,  a  radius  of  0.25  m  and  its 
centre  is  located  1  m  below  the  water/air  interface.  For  the  cylinder,  there  will  be  a 
combination  of  specular  rehections  from  the  length  of  the  cylinder,  endcap  rehections, 
and  diffractions  from  the  edges  of  the  endcaps.  In  Fig.  4,  we  show  the  computed 
spectra  for  the  3  angles  of  source/receivers  (the  same  angles  as  for  the  sphere)  for  all 
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interactions  (solid)  and  the  single  scatter  solution  (markers).  As  with  the  sphere,  the 
spectra  for  the  60°  incidence  are  very  similar,  for  30°  there  are  observable  differences, 
and  for  normal  incidence  the  differences  are  large. 

It  is  interesting  to  consider  the  solution  of  the  BIEM  without  using  the  additional 
interior  points  for  Schenke’s  method.  From  Eq.(18),  the  values  of  the  irregular  fre¬ 
quencies  can  be  predicted  for  the  case  of  flat-ended  cylinder.  In  Fig.  5a,  we  show 
the  computed  spectra  (a  5  Hz  spacing  is  used)  for  the  3  source/receiver  angles  used 
in  the  previous  example.  Both  the  azimuthally-uncoupled  solution  (single-scatter, 
dashed)  and  the  fully  coupled  solution  (solid)  are  shown.  However,  no  additional 
interior  points  are  used  in  the  computations.  As  can  be  seen  there  are  frequencies 
near  which  the  spectral  amplitude  is  poorly  behaved  as  a  function  of  frequency.  This 
behaviour  is  more  noticeable  for  some  source/receiver  angles.  This  is  likely  due  to 
the  differing  projections  of  the  solution  onto  the  eigenfunction  corresponding  to  this 
frequency.  The  frequency  values  predicted  by  Eq.(18)  are  shown  in  Fig.  5a  as  vertical 
lines.  In  Fig.  5b,  the  computed  Condition  Number  (ratio  of  amplitude  of  largest  to 
smallest  eigenvalue  for  a  matrix  using  the  MATLAB  routine)  for  the  fully  coupled  ma- 
trix(solid)  and  uncoupled  matrix  (dashed,  single-scatter  approximation)  are  shown. 
The  condition  number  behaviour  is  so  similar  for  the  coupled  and  uncoupled  matrices 
that  is  difficult  to  distinguish  the  different  lines  in  Fig.  5b.  Once  again,  for  the  full 
and  single-scatter  method,  the  poor  solution  behaviour  is  related  to  the  proximity  of 
irregular  frequencies  and  as  discussed  earlier,  these  values  are  the  same  for  the  fully 
coupled  and  uncoupled  matrices. 

Finally,  we  consider  the  case  of  the  cylinder  tilted  at  -45°  (see  Fig.  2)  with  respect 
to  the  vertical.  Three  different  source/receiver  positions  are  considered,  10  m  range 
and  at  45°,  60°,  and  90°  off  the  vertical  (and  to  the  right  of  the  tilted  cylinder  in 
Fig.  2)  .  For  the  source/receiver  at  45°  there  is  a  specular  reflection  off  the  endcap  of 
the  cylinder.  For  the  source/receiver  at  60°  off  vertical,  there  is  a  ray  which  reflects  off 
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Figure  4:  A  comparison  of  the  computed  backscattered  spectra  using  the  cylindrical 
BIEM  with  the  full  solution  (solid  lines)  and  the  single-scatter  solution  (lines  with 
markers)  for  source/receiver  60°  (circles),  30°  (diamonds),  and  0°  (+)  off  the  vertical 
axis  below  the  cylinder.  For  clarity,  the  symbols  are  shown  at  every  fifth  frequency 
value. 
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Figure  5:  (a)  The  computed  backscattered  spectral  amplitude  for  the  source/receiver 
configurations  of  Fig.  4  using  full  solution  (solid)  and  single-scatter  solutions  (dashed). 
There  have  been  no  additional  interior  points  used  in  the  computations,  (b)  The 
computed  matrix  condition  number  for  the  full  solution  matrix  (solid)  and  uncoupled 
single-scatter  matrix  (dashed).  For  both  plots,  the  vertical  lines  indicate  the  predicted 
locations  of  the  eigenfrequencies  of  the  interior  Dirichlet  problem. 
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the  upper  surface  and  is  then  incident  upon  the  cylinder  at  close  to  broadside.  This 
energy  is  then  reflected  back  off  the  top  surface  before  returning  to  the  receiver.  This 
is  still  a  single-scattering  event  as  there  is  only  one  reflection  off  the  cylinder.  Finally, 
for  the  source/receiver  along  the  horizontal  (90°),  the  tilted  cylinder  is  a  type  of 
corner  reflector.  A  horizontal  ray  is  reflected  vertically  upwards  and  then  is  reflected 
by  the  upper  surface  downwards  and  is  then  converted  back  to  a  horizontal  ray  after 
a  second  reflection  from  the  cylinder.  In  Fig.  6  we  show  the  computed  spectra  with 
all  interactions  (solid  line)  and  with  the  single  scatter  approximation  (markers).  The 
circles  correspond  to  the  45°  source/receiver,  the  diamonds  the  60°  case,  and  the 
plus  signs,  the  horizontal  source/receiver.  The  largest  differences  between  the  single¬ 
scatter  (line  with  +  markers  at  the  bottom  of  the  figure)  and  the  full-scattering 
solution  (unassociated  line  near  the  middle  of  the  figure)  occur  for  the  horizontal 
source/receiver.  This  large  difference  is  anticipated  from  the  previous  discussion  with 
respect  to  the  strong  echo  from  2  interactions  with  the  cylinder.  In  Fig.  7  we  show  the 
computed  echo  time  series  for  an  incident  2500  Hz  Ricker  wavelet.  The  large  multiple- 
scattering  echo  can  be  seen  in  Fig.  7c  as  predicted  by  simple  ray  considerations  (solid 
lines,  full  solution,  dashed  lines,  single-scatter  solutions). 
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Figure  6:  A  comparison  of  the  computed  backscattered  spectra  using  the  cylindrical 
BIEM  with  the  full  solution  (solid  lines)  and  the  single-scatter  solution  (markers)  for 
source/receiver  45°  (circles),  60°  (diamonds),  and  90°  (+)  off  the  vertical  axis  below 
the  cylinder.  The  cylinder  is  tilted  at  -45°  with  respect  to  the  horizontal.  For  clarity, 
the  symbols  are  shown  at  every  fifth  frequency  value.  The  full-scattering  solution  for 
90°  is  the  line  near  the  vertical  centre  of  the  plot  which  appears  unassociated. 
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Figure  7:  The  computed  echo  time  series  for  a  2500-Hz  Ricker  wavelet  using  the 
full  solution  (solid  line)  and  the  single-scatter  solution  (dashed)  for:  (a)  45°  off  the 
vertical,  (b)  60°  off  the  vertical  and  (c)  90°  off  the  vertical. 
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IV.  DISCUSSION  OF  RESULTS 

In  this  paper  a  modified  azimuthally-symmetric  BIEM  was  presented  for  the  solu¬ 
tion  of  scattering  from  cylindrical  objects  below  a  pressure-release  interface.  Although 
the  results  shown  are  for  a  sphere  and  a  flat  end-capped  cylinder,  the  methods  de¬ 
scribed  are  valid  for  more  general  target  shapes  with  azimuthal  symmetry.  The  scat¬ 
tering  process  for  the  cylindrical  object  in  free  space  can  be  constructed  by  solving 
a  sequence  of  two-dimensional  problems,  one  for  each  angular  order.  When  the  in¬ 
teractions  between  the  object  and  the  interface  are  considered  then  there  is  coupling 
between  the  Fourier  orders.  The  computation  of  this  azimuthal  conversion  matrix  is 
described  for  the  cases  of  the  cylinder  lying  parallel  to  the  interface  and  also  for  the 
cylinder  tilted  with  respect  to  the  interface.  For  the  case  of  the  parallel  cylinder,  this 
conversion  matrix  can  be  efficiently  computed  by  using  a  wavenumber  integral.  When 
the  conversion  matrix  is  set  to  zero,  the  resulting  computed  pressure  field  corresponds 
to  the  single-scattering  approximation. 

In  the  numerical  examples,  the  proposed  method  was  first  benchmarked  against 
another  solution  method  for  a  sphere  and  the  agreement  was  found  to  be  excellent. 
Different  scattering  geometries  were  considered  for  a  cylinder  parallel  to  the  interface 
and  for  a  cylinder  tilted  at  -45°  with  respect  to  the  vertical.  For  the  parallel  cylinder 
it  was  found  that  for  small  monostatic  scattering  angles  with  respect  to  the  cylinder’s 
axis,  the  single-scattering  solution  was  accurate  and  became  less  accurate  as  the 
scattering  angle  increased.  For  the  tilted  cylinder,  the  cylinder  and  the  interface 
combined  to  form  a  corner-like  reflector  and  in  this  case  there  was  a  strong,  multiple¬ 
scattering  return  for  the  source/receiver  parallel  to  the  interface.  The  significant 
difference  with  the  single-scattering  solution  was  evident  in  both  the  spectral  and 
time  domains. 

In  the  method  of  this  paper,  the  problem  of  irregular  frequencies  was  addressed  by 
using  additional  interior  equations  in  the  BIEM.  In  order  to  investigate  the  location 
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of  the  irregular  frequencies  for  the  single-scatter  (uncoupled)  and  fully-coupled  ap¬ 
proaches  the  extra  interior  equations  were  removed.  The  computational  values  of  the 
irregular  frequencies  appeared  to  be  in  excellent  agreement  with  the  theoretically  pre¬ 
dicted  values  and  these  values  were  the  same  for  the  single-scatter  and  fully-coupled 
formulations. 
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List  of  Figures 

Figure  1  -  Schematic  drawing  of  a  1.5m  long  cylinder  with  radius  0.5.  A 
generic  receiver  position  is  shown  at  ( r,6,z  =  0).  A  point  on  the  surface  at 
(r'  =  0.5, 9',  z'  =  0)  is  also  shown. 

Figure  2  -  Schematic  drawings  of  geometry  and  notation  for  (upper)  tilted  cylin¬ 
der  with  a  source  and  the  mirror  image.  The  tilted  cylinder  has  a  centre  offset 
H.  (lower)  The  angular  coordinate  system  is  shown  on  the  endcaps  of  a  cylinder 
lying  parallel  to  the  interface.  A  receiver  point  P  is  indicated  with  the  angles 
6r  and  9$  the  azimuthal  angles  measured  from  the  lower  and  image  cylinders, 
respectively. 

Figure  3  -  A  comparison  of  the  computed  backscattered  spectra  using  the 
virtual  source  technique  (sold  lines)  and  the  cylindrical  BIEM  (markers)  for 
source/receiver  60°  (circles),  30°  (diamonds),  and  0°  (+)  off  the  vertical  axis 
below  the  sphere.  For  clarity,  the  symbols  are  shown  at  every  fifth  frequency 
value. 

Figure  4  -  A  comparison  of  the  computed  backscattered  spectra  using  the  cylin¬ 
drical  BIEM  with  the  full  solution  (solid  lines)  and  the  single-scatter  solution 
(lines  with  markers)  for  source/receiver  60°  (circles),  30°  (diamonds),  and  0°  (+) 
off  the  vertical  axis  below  the  cylinder.  For  clarity,  the  symbols  are  shown  at 
every  fifth  frequency  value. 

Figure  5  -  (a)  The  computed  backscattered  spectral  amplitude  for  the 
source/receiver  configurations  of  Fig.  4  using  full  solution  (solid)  and  single¬ 
scatter  solutions  (dashed).  There  have  been  no  additional  interior  points  used 
in  the  computations,  (b)  The  computed  matrix  condition  number  for  the  full 
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solution  matrix  (solid)  and  uncoupled  single-scatter  matrix  (dashed).  For  both 
plots,  the  vertical  lines  indicate  the  predicted  locations  of  the  eigenfrequencies 
of  the  interior  Dirichlet  problem. 

Figure  6  -  A  comparison  of  the  computed  backscattered  spectra  using  the  cylin¬ 
drical  BIEM  with  the  full  solution  (solid  lines)  and  the  single-scatter  solution 
(markers)  for  source/receiver  45°  (circles),  60°  (diamonds),  and  90°  (+)  off  the 
vertical  axis  below  the  cylinder.  The  cylinder  is  tilted  at  -45°  with  respect  to 
the  horizontal.  For  clarity,  the  symbols  are  shown  at  every  fifth  frequency  value. 
The  full-scattering  solution  for  90°  is  the  line  near  the  vertical  centre  of  the  plot 
which  appears  unassociated. 

Figure  7  -  The  computed  echo  time  series  for  a  2500-Hz  Ricker  wavelet  using  the 
full  solution  (solid  line)  and  the  single-scatter  solution  (dashed)  for:  (a)  45°  off 
the  vertical,  (b)  60°  off  the  vertical  and  (c)  90°  off  the  vertical. 
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Figure  1:  Schematic  drawing  of  a  1.5m  long  cylinder  with  radius  0.5.  A  generic 
receiver  position  is  shown  at  (r,  9,  z  =  0).  A  point  on  the  surface  at  (rf  =  0.5,  9',  z'  =  0) 
is  also  shown. 
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Figure  2:  Schematic  drawings  of  geometry  and  notation  for  (upper)  tilted  cylinder 
with  a  source  and  the  mirror  image.  The  tilted  cylinder  has  a  centre  offset  H.  (lower) 
The  angular  coordinate  system  is  shown  on  the  endcaps  of  a  cylinder  lying  parallel  to 
the  interface.  A  receiver  point  P  is  indicated  with  the  angles  9r  and  9,  the  azimuthal 
angles  measured  from  the  lower  and  image  cylinders,  respectively. 
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Figure  3:  A  comparison  of  the  computed  backscattered  spectra  using  the  virtual 
source  technique  (sold  lines)  and  the  cylindrical  BIEM  (markers)  for  source/receiver 
60°  (circles),  30°  (diamonds),  and  0°  (+)  off  the  vertical  axis  below  the  sphere.  For 
clarity,  the  symbols  are  shown  at  every  fifth  frequency  value. 
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Figure  4:  A  comparison  of  the  computed  backscattered  spectra  using  the  cylindrical 
BIEM  with  the  full  solution  (solid  lines)  and  the  single-scatter  solution  (lines  with 
markers)  for  source/receiver  60°  (circles),  30°  (diamonds),  and  0°  (+)  off  the  vertical 
axis  below  the  cylinder.  For  clarity,  the  symbols  are  shown  at  every  fifth  frequency 
value. 
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Figure  5:  (a)  The  computed  backscattered  spectral  amplitude  for  the  source/receiver 
configurations  of  Fig.  4  using  full  solution  (solid)  and  single-scatter  solutions  (dashed). 
There  have  been  no  additional  interior  points  used  in  the  computations,  (b)  The 
computed  matrix  condition  number  for  the  full  solution  matrix  (solid)  and  uncoupled 
single-scatter  matrix  (dashed).  For  both  plots,  the  vertical  lines  indicate  the  predicted 
locations  of  the  eigenfrequencies  of  the  interior  Dirichlet  problem. 
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Figure  6:  A  comparison  of  the  computed  backscattered  spectra  using  the  cylindrical 
BIEM  with  the  full  solution  (solid  lines)  and  the  single-scatter  solution  (markers)  for 
source/receiver  45°  (circles),  60°  (diamonds),  and  90°  (+)  off  the  vertical  axis  below 
the  cylinder.  The  cylinder  is  tilted  at  -45°  with  respect  to  the  horizontal.  For  clarity, 
the  symbols  are  shown  at  every  fifth  frequency  value.  The  full-scattering  solution  for 
90°  is  the  line  near  the  vertical  centre  of  the  plot  which  appears  unassociated. 
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Figure  7:  The  computed  echo  time  series  for  a  2500-Hz  Ricker  wavelet  using  the 
full  solution  (solid  line)  and  the  single-scatter  solution  (dashed)  for:  (a)  45°  off  the 
vertical,  (b)  60°  off  the  vertical  and  (c)  90°  off  the  vertical. 
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